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I . Introduction 


In  this  paper  we  consider  approximation  techniques  for 
functional  differential  equations  (FDE's)  based  on  classical 
"least  squares"  or  best  spline  approximations.  In  a recent 

paper  [ 7 ] , Burns  and  Cliff  developed  an  approximation  scheme 
that  employed  piecewise  linear  function  approximations.  However, 
to  our  knowledge,  the  ideas  for  use  of  splines  (piecewise  linear 
or  higher  order)  in  approximating  FDE's  as  developed  in  our  pre- 
sentation here  are  new.  We  develop  our  ideas  in  the  context  of 
an  abstract  approximation  theorem  (the  Trotter-Kato  approximation 
theorem  in  linear  semigroup  theory)  which  greatly  facilitates 
arguments  to  establish  convergence.  Use  of  the  Trotter-Kato 
theorem  in  this  way  is  not  new;  indeed,  it  has  been  used  recently 
in  connection  with  other  approximation  schemes  for  problems  in- 
volving linear  and  nonlinear  FDE's  (see  [3  ] for  a survey  and  in 
addition,  the  recent  papers  [1], [2] , [4],  [7],  [8],  [9]). 

We  restrict  our  considerations  in  this  paper  to  linear 
systems  and  their  approximation.  The  use  of  the  ideas  presented 
here  in  nonlinear  system  problems  along  with  applications  to  op- 
timal control  and  parameter  estimation  problems  will  be  discussed 
elsewhere.  Our  main  purpose  here  is  to  develop  the  fundamental 
theoretical  ideas  for  spline  based  methods  and  present  a sample 
of  our  findings  in  numerical  experiments  with  these  approximations. 
As  will  be  evident  from  our  discussions  of  the  numerical  results 
in  Section  5 below,  we  believe  that  these  spline  based  approxima- 
tion techniques  can  offer  significant  advantages  over  other  methods 


A 


2. 


PM 

tl 


(e.g..  those  discussed  in  [1],  [3],  [7],  [8]  and  [9])  in 
many  instances. 

We  first  develop  in  Section  2 a general  setting  for  our 
system  approximation  problem  in  an  appropriate  Hilbert  space.  ! 

Much  of  the  material  in  this  section  is  closely  related  to  known 
results,  but  must  be  presented  in  order  to  give  a complete  and 
clear  discussion  of  our  ideas.  A general  approximation  result 
is  developed  in  Section  3 and  we  then  in  Section  4 show  how  spline 
based  methods  can  be  treated  in  a simple  manner  as  a special  case 
of  the  approximation  ideas  of  Section  3. 

One  feature  of  our  presentation  is  that  the  proofs  in 
Section  4 follow  immediately  from  standard  results  in  spline  theory. 

Aesthetically,  our  development  has  appeal  since  the  theoretical 
foundations  are  cleaner  than  those  for  other  types  of  methods  pro- 
posed to  date  (e.g.  see  [3],  [7]).  From  a practical  viewpoint 
these  methods,  as  we  have  already  noted,  should  prove  advantageous 
in  many  control  and  identification  problems. 

The  following  notation  will  be  used  throughout  the  paper. 

For  -»  < a < b < °°,  L2(a,b;Rn)  is  the  Hilbert  space  of  equivalence 
classes  of  all  functions  x:[a,b]  -*■  Rn  such  that  | x | ^ is  in- 
tegrable.  |*|  is  the  Euclidean  norm  in  Rn,  (•,•) 2 and  1‘M 
denote  the  usual  inner  product  and  norm  in  I^fa.bjR11).  C^(a,b;Rn), 
k = 0,1,2 denotes  the  space  of  Rn  -valued  continuous  func- 


tions which  possess  k continuous  derivatives  on  [a,b].  (At  end 
points  of  closed  intervals,  "derivative"  of  course  means  the  ap- 


3. 


propriate  one  sided  derivative.)  For  k = 0 this  is  the  usual 

space  of  continuous  functions  which  wc  denote  simply  by  C(a,b;Rn). 

1 7 n 

W ’ (a,b;R  ) is  the  space  of  absolutely  continuous  functions 
x:[a,b]  -*•  Rn  such  that  x G I^fa.bjR0).  The  space  of  equivalence 
classes  of  all  functions  x:[0,°°)  -*■  Rn  such  that  | x | ^ is  in- 
tegrable  on  bounded  intervals  will  be  denoted  by  1^  joc  = 

L 2 }oc  (0  »°°  ; Rn)  . For  any  t^  > 0 the  collection  of  all 
f G L-,  restricted  to  [0,t^j  gives  the  space  L^CO.t^R11), 

of  course.  In  the  special  case  where  a = - r , b = 0 with  r > 0 
we  shall  abreviate  the  notation  introduced  above  and  simply  write 
L^,  C , etc.  The  state  space  for  our  considerations  will  be  the 
Hilbert  space  Z = Rn  * with  norm  | Cn » d> ) l ^ = (|n|2  + 

(n,4>)  6 Z , and  inner  product  < (n1 , ) , (n2  » >Z  = nln2  + < *1  »*2  ^ ? ’ 

C n i * ^ i ) ^ i = 1,2.  Corresponding  to  C , k = 0,1,2,...,  we 
introduce  the  linear  subspaces  of  Z defined  by 

if  - { (4>  (0)  , 40  | 4>  CC^}.  Again,  for  k = 0 we  write  if  instead 
of  Throughout  this  paper  we  shall  not  distinguish  notationallv 

between  a function  4>  and  its  equivalence  class  in  L^.  If  4>  is 

A 

in  C,  then  4>  denotes  the  element  (4>(0),4O  in  if.  If  on  the 

A 

other  hand  <)>  is  an  element  in  if,  then  4>  denotes  the  unique 

A 

element  in  C such  that  4>  = (4>  (0)  , 4> ) . The  identity  map  in  any 
space  and  its  matrix  representations  will  always  be  denoted  by  I. 

For  a function  x:[-r,a)  -*•  Rn,  a > 0,  the  symbol  xt, 
tG[0,a),  denotes  the  function  [-r,0]  -*•  Rn  defined  by 
xt(0)  = x(t+0),  0G[-r,O).  Finally,  for  any  function  x 


of  one 


interchangeably  e i t h e r 


x or 


4 . 


independent  variable  we  < !i a 1 l use 
!'\  to  denote  the  derivative  of  \ 


with  respect  to  this  variable. 


2.  rite  linear  1 unc  t i ona  1 differential  equation 

l it  this  section  we  detail  the  type  of  equation  we  consider 
in  this  paper  and  state  the  results  for  this  equation  which  are 
important  for  the  developments  in  the  following  sections. 

(>iven  ( 2 and  t (.  I.,  j we  consider  the  Cauchy 

p rob  l em 


x ( t ) - 1. 1 x 1 1 ♦ f ( t ) , t ^ 0, 
( x ( 0 ) , x Q ) - ( n , $ ) . 


(2.1) 


lor  simplicity  of  exposition  we  assume  that  L (x 1 1 lias  the  form 

m ,0 

l.  (x  ) - l A.x(t-T.)  ♦ A ( 0 ) x ( t ♦ 0 ) d 0 (2.2) 

1 i » 0 1 1 i - r 

where  0 - t ^ <...<  r and  A.,A(0)  are  n * n matrices, 

the  elements  of  0 * M0)  being  square  integrahle  on  ( r,0].  Not 
that  tor  anv  right  hand  side  of  tvpe  (2.2)  we  can  assume  t = r 
(possible  after  redefining  A(0)  or  with  A ■ 0).  The  form  of 
I ( x t ) given  in  (2.2)  is,  to  our  knowledge,  sufficiently  general 
to  include  all  linear  autonomous  1 PI  ' s arising  in  applications. 

(In  regard  to  more  general  equations  see  the  remark  at  the  end 
of  this  sec t ion . ) 


5. 


A solution  of  (2.1)  is  a function  x:[  r,u)  •*  Rn , a > 0, 
which  is  absolutely  continuous  on  [ 0 , u ) such  that  (x(.0),Xq)  = 
(n,f)  and  for  0 < t s a, 


x (t) 


0 


L(.\  )do  f 


f (o ) do , 


C.3) 


A solution  of  (2.1)  is  denoted  bv  x(t)  x(t;n,<f,f).  In  case 
f 0 we  shall  write  simply  x(t;n,<f). 


Note,  that  for  any  two  functions  x,v:[-r,o.] 


. u j > 0 , 


which  are  absolutely  continuous  on  and  satisfy  x = y^ 

in  L7  and  x(t)  = y(t)  for  t fc  [0,a.],  the  maps  o L(x  ) 
and  o -*■  l.(y  ) are  defined  and  in  the  same  equivalence  class  of 
Rn-valued  functions  on  [0,o.|.  In  general  1 is  not  defined  on 
all  of  Z,  but  can  be  defined  on  r-  by  L(4>)  = 1(f),  if  c C. 

The  following  lemma  is  stated  without  its  quite  elementary 

proo f . 

I.emma  2.1.  Ihe  solutions  of  (2.1)  exist  on  ( r,v)  and  are  unique 

lv  determined.  Moreover,  for  sequences  ♦ (n,f)  in  2 

and  f.  -►  f in  ! , , we  have  for  all  T s 0 

k — 2 , loc  


sup  | x ( t ) - x.  ( t ) | -*•  0 
t t ( 0 ,T ] k 


as  k •*  ™ where  x(t)  = x(t;n,$,f)  and  x^(t)  = x(-t:nk,'}'k,tk^‘ 


In  the  case  of  the  homogeneous  equation,  i.e. 


f 0 , w e d e f i n e 


6. 


the  family  S(t),t  > 0,  of  operators  Z - Z by 


S ( t ) Cn , 4> ) = (x(t;n,«J>),xt(n,*))  (2-4) 


for  each  (n,$)  £ 2. 

Lemma  2.2. 

a)  S(t),t  > 0,  is  a C - semigroup  of  bounded  linear 
operators . 

b)  The  infinitesimal  generator  jy  of  this  semigroup 
and  its  domain  &{sst ')  are  given  by 

i?(V)  = ( (<t>  (0)  ,*)  G Z|*  6 W1*2}, 
j y(<p(0),<p)  = (L(4>),4>),  (4>(o)  ,<j>)  g • 

c)  For  k = 1,2,...  the  sets 

= U Gi5k|K0)  = L(4») } 

and 

for  A sufficiently  large 

are  dense  in  Z . 

Parts  a)  and  b)  of  this  lemma  are  by  now  standard  results 
if  one  deals  with  FDF's  in  the  state  space  Z (cf.  for  instance 
[3],  [91,  [15]). 


we  first 


In  order  to  establish  density  of  ( wv  - V l ) — we  first 
observe  that  t>  , k = 0,1,2,...,  is  dense  in  Z because 

a k)  C r- K 1 and  £/(  vy  ) = Z for  k = 1,2,...  . Ibis  last 

fact  is  true  for  the  infinitesimal  generator  of  any  C -semigroup 
of  bounded  linear  operators  in  Z (see  for  instance  any  standard 


reference  on  linear  semigroup  theory  such  as 


or  ( 1 0 ] ) . Then 


density  of  (a  \ 1 ) _ * ^ for  \ sufficiently  large  follows  once 
one  argues 


( ^ -M)-' 


(2.5) 


for  k = 1,2,...  and  such  X. 

The  resolvent  operator  (a  XI)  ^ exists  for  X suf- 
ficiently large  and  is  a bounded  linear  operator  Z ^_’(a). 
Therefore,  given  r G r- ^ the  equation 


( . a - A I ) 4>  = v 


(2.0) 


has  a unique  solution  41  for  X sufficiently  large.  We  only 

A g 

have  to  check  it  41  <=  L/  . but  using  the  definition  of  ci  we 
see  that  (2.0)  is  equivalent  to 


4 - X4  = v and  L(4)  - X4>(0)  = f(0) 


which  implies  4>  G CK  and  MO)  = \4>(0)  ♦ V(0)  = L($).  This 


8. 


proves  (2.5).  Finally,  density  of  1/^  now  follows  at  once  from 
_/  = ( .i  -XI)  * r ^ * for  X sufficiently  large. 

A fundamental  notion  we  shall  use  in  our  approach  is  the 
concept  of  dissipativeness.  Recall,  that  a linear  operator 
B:  _/(B)  -*■  Z,  C Z,  is  called  dissipative  in  Z if  (Bz,z'-  < 0 

tor  all  z t £/(B).  If  B is  the  infinitesimal  generator  of  a 
f ■ semi  group  T(t),  t ^ 0,  of  bounded  linear  operators  Z Z and 
H - to  I is  dissipative  for  some  io  f R,  then  | T ( t ) z | - £ ewt  I z |, 

tor  all  : ( Z and  t _>  0.  The  solution  semigroup  S(t),  t > 0, 

detined  in  (2.4)  in  general  satisfies  an  estimate 

I S (t)  z | z < Me(JJt  | z | : , z € Z,  t > 0, 

where  to  € R and  M > 1.  Therefore  a - col  will  not  he  dissipa- 

tive in  Z for  any  to  £ R.  Following  an  idea  used  in  [15)  we 
introduce  an  equivalent  norm  in  Z such  that  Z with  this  norm 

is  again  a Hilbert  space  Z and  for  some  constant  w E R 

S 

y 10 1 is  dissipative  in  the  space  Z . Roughly  speaking, 
the  new  norm  reflects  the  weight  placed  by  (2.2)  on  different 
parts  of  the  past  history. 

Corresponding  to  the  difference  part  in  (2.2)  we  define 


the  weighting  function  g to  be  a step  function  on  [ - r , 0 ) 
such  that 


The  space 


supplied  with  the  norm 


| (n  ,4>)  | 7 = ( | n | 2 - 1 4>(oj  |2g(0)de)1/z,  (n,4>)  e Z, 


will  be  denoted  by  Z . Obviouslv,  Z is  a Hilbert  space  with 

g g 

inner  product 


( (n  i ,$i ) , (n->  ) >7  = rh  n 7 


$,(0)g(d)dd, 


t z,  i = 1 ,2 


Following  the  arguments  given  in  [3,  p.  186],  we  can  establish  the 


following  result. 


Lemma  2.3.  A - wl  is  dissipative  in  Zn,  ie, 

o 


(Si/z  , z )z  <_  to  i z | ^ 
g g 


for  z f7(.y), 


wi  th 


m+  1 

to  = 7--  + 


'V  * 1 jjhl2  * 7 fjMe)! 


Obviously,  the  norm  |*|7  is  equivalent  to  | - | ^ • 

g 

Therefore,  all  results  stated  in  Z which  involve  only  topological 

concepts  remain  valid  in  Zi  (e.g.,  the  results  of  Lemma  2.2). 

g 

Lemma  2.3  will  be  basic  for  our  use  of  the  Trotter-Kato 


Theorem  in  the  next  section  in  order  to  obtain  an  approximation  result 


10. 


for  the  homogeneous  Cauchy  problem,  i.e.  (2.1)  with  f = 0. 

We  now  turn  to  the  general  nonhoinogeneous  case  of  (2.1) 
and  define 

■t 

z(t)  = z(t;n,<t>,f)  = s (t ) (n  , 4>)  + S (t-o)  (f  (a)  , 0)  do  , t > 0 (2.7) 

J0 

for  (n,<J>)  6 Z and  f FL^  loc.  In  applications  (especially  con- 
trol problems)  it  is  also' important  to  investigate,  for  fixed 
tx  > 0,  the  operator  (0 , t ; Rn)  - C(0,t]L;Z)  defined  for 

0 < t < tj  and  f G L2(0,t1;Rn)  by 

(^f)(t)  = J*  S(t-o)  (f(o)  , 0 ) do.  (2.8) 

lor  our  discussion  below,  the  first  part  of  the  following  lemma  is 
an  essential  equivalence  result.  The  second  part,  while  not  es- 
sential to  our  presentation  here,  is  needed  in  applying  our  approxi- 
mation results  to  control  problems  (see  [ 2]),  and  is  stated  for 
future  reference. 

Lemma  2.4. 

a ) Civen  ( n , <f> ) £ Z and  f G L 2 j 0 c » we  have 
z(t)  = (x(t;n,<f>,f)  ,xt(n,4>,f))  for  t^0 


1 


u a c i c 


*.  l.  l ) tb  uoiinea  o\  ( i . / j . 


b ) for  any  t ^ > 0 the  operator  ■/  defined  by  (2.8) 
is a compact  1 ijiea r operator . 


The  results  of  this  Lemma  have  already  been  established 
elsewhere  (see  Theorems  2.1  and  3.2  of  [ 2 ],  Theorem  1 of  [4], 

I.emma  1.3.7  of  [11]  or  Lemma  2.2  of  [ 12])  but  for  the  sake  of  com- 
pleteness we  shall  sketch  £he  arguments  needed  to  give  a proof. 

\ , 

We  first  assume  (^0),<fc)  EMa/)  and  f GC^O.-iR"). 

Then  (see  for  instance  [5  , p.  31])  z(t)  given  in  (2.7)  is  the 
unique  strong  solution  to  the  abstract  Cauchy  problem 


z(t)  = tvz(t)  + (f(t)  ,0)  , t > 0, 

2(0)  = U(0),<t>). 


(2.9) 


On  the  other  hand  for  x(t)  = x (t ; <f>  (0)  , 4> , f ) we  have 

x 6 1^’  (-r,tpRn)  for  any  tj  > 0,  which  implies  x^.  E W*,“ 

for  all  t _>  0.  By  a characterization  of  functions  in 

(see  [S,  p.  2 1 | ) , the  derivative  of  the  map  t -+  x( 

d 

considered  as  a map  [ 0 , 00 ) -*■  L,  exists  and  t = where 

x (0)  = x(t  + 0),  0E[-r,0].  This  together  with  (2.1)  reveals  that 
(x(t),xt)  is  also  a solution  of  (2.9).  By  uniqueness  of  solutions 
of  (2.9)  we  have  z(t)  = (x(t),x.),  t 0 and  thus  the  equiva- 
lence in  part  a)  of  Lemma  2.4  obtains  whenever  (n,4>)  £ £^(jy) 


• *••»**  ‘ * — i - ■ • ' ' • 'S  ■ «*r 


--  •- 


\ 


,1 


and  f GC*(0,»;Rn). 


To  complete  the  proof  of  part  a)  we  observe  that 


y ( \y~)  - Z and  C (0,t^;Rn)  is  dense  in  L-,(0,t^;Rn)  for  any 
t^  >0.  Moreover,  for  0 <t£t ^ , z(t)  given  by  (2.7)  and  (x(t),x) 
x(t)  = x (t ; n , <J> , f ) , depend  continuously  on  (n,4>)  G Z and 
f € L2(0,ti:Rn). 

In  order  to  prove  part  b)  we  fix  t ^ 0 and  a bounded 

subset  (I  of  L-,  (0  , t Rn)  . We  have  to  prove  that,  for  each  t 
in  [O.tjJ,  the  set  {(  rf)(t)|f  6G)  is  precompact  in  Z and 
that  ( '*'f|f  G G}  is  bounded  and  equicontinuous . 

If  we  observe  that  ( y f) (t)  = (x(t;0,0,f),  xt(0,0,f)), 
tGlO.tj],  and  that  {x  ( • ; 0 , 0 , f ) | f € G}  and  {x  ( • ; 0 ,0  , f)  | f G GJ 
arc  bounded  subsets  of  C(-r,tjjRn)  and  L , (- r , t . ; Rn)  , res- 
pectively, it  is  immediate  that  { ( Vf ) (t)  | f G G)  is  precompact 

in  Y,'  endowed  with  the  norm  j <#> ) = sup  | 4>  | and  therefore  also 

# l-r ,0] 

in  Z.  The  same  arguments  also  give  equicontinuity  and  boundedness 
of  { "/f  | f G G) . 


Remark . All  rcsuLts  given  in  this  section  remain  valid  for  more 
general  equations  than  considered  here.  L must  only  be  a continuous 
functional  from  G to  Rn  satisfying  the  conditions  given  by 

Bor i sov ic -Turbab i n in  [6]  (sec  also  [2]).  If 

f° 

1.(40  = [dp  (0 ) ] <)>  (0)  , (j>  G C,  then  the  weighting  function  g in 

* - r rd 

this  more  general  case  has  to  be  defined  by  g(0)  = 1 + | dn | , 

0G[-r,0], 


} ■ A general  approximation  scheme . 

Masic  to  our  approach  is  the  following  version  of  the 
so-called  Trotter-Kato  Theorem  (see  | 1 (.) , Thm.  4.t>}). 


Lem  in  a 3.1.  Let  T(t)  and 

semigroups  in  a Banach  space 
J and  eg  , re  spec  t i ve  1 y . 


, N C t ) , N = 1,2, . . . , t > 0 be  C - 

— — o 

A infinitesi ma  1 gene ra tors 

Assume  that  the  lol lowing  condi t ions 


are  satisfied: 


i)  (§ t ab i lit  y hypo  the s is).  There  exist  constants 

w HOl  WN  such  that  g-  UI  and  ^ N - ^ [ are 

dissipative  on  V and  the  sequence  f coN } is  hounded. 

(Lons i s t ency  hypo  thesis).  There  exists  a subset 

^ Ls  ( ft ) O £/(  ft  ) which  together  with 

N-l  

( WAI)^  for  some  A - 0 is  dense  in  V and  such 

that  ft  -►  gy  Co r all  y G £A 

Then 

lim  TN(t)y  = T(t)y 

N->a> 


for  jil  l 


y A uiiiformlv  on  hounded  t- intervals. 


We  call  ( Z ,1’  , W },  N = 1,2 an  approximation 

scheme  for  the  Cauchy  problem  (2.1)  if  {ZN}  is  a sequence  of 
subspaces  ot  2^,  {l’N}  is  the  sequence  of  orthogonal  projections 

*’  :Zg  *■  2 an(l  (tv'1^}  is  a sequence  of  operators  Z -*•  ZN. 


••-a  • > . ; .• 


v.< 


!«.»>•  »»  • >*.  . 


"'.S  . 


*•  ...  I 


Theorem  3.1. 

Suppose 

that 

(ZN,p\yN) 

is  an  approximation 

scheme  for  ( 

2 . 1 1 satis 

tying 

t he  f o 1 1 o\»  i ng 

cond i t i ons : 

U)  ZN 

C y) , 

N = 

1 2 

*»*’»•••  • 

li  i ) lv  N 

■ pVp", 

N = 

1 2 

*»•■»•••  • 

(iii)  (a) 

1 im  l'N; 

N 

; ; in  2 for  a 1 1 

z t Z , 

(b) 

Lor  some 

into 

g e r k > I u o 

have  1 im  l.  ( ig‘  ) = I.  (*’ ) 

„n 

i n R 

where  y 

and 

N 

1 s 

1 im  IV  IV 

N *•» 

defined  by  P 

N .... 

in  1 , tor  all  V € (' 

N ' N N 

* (V  (O')  . V ^ 1 . 

Then 

each  y 

is  the  infinitesimal  generator  of  a 

t)  - sem  i g roup 

sN(t),  r 

" 0, 

such  that 

s 

V )Z 

N _N  M 

C - , N 

1,2 

a ml 


1 im  SN(t)z  = S(t):  (.3.1') 

N ►“» 


for  a 1 1 z C Z,  uniformly  on  bounded  t intervals. 


Proof . We 

in  Z wit h 
S 

map  coincide 
ZN  C £/(y)  , 


first  show  that  or N - wl , N = 1,2,,..,  is  dissipative 

N 

to  as  given  in  Lemma  2.3.  Since  1'  and  its  dual 

N 

, we  get,  by  using  the  definition  of  y and 
the  following  estimate 


M . -I*.  t \ -A-**  * 


**“rry  rtvrv-  •vjk'**.  . a*  4.  . 


• « .m 


■4 


*!•»  ' * HIsav  Jft‘. 


for  all 


(LyNz,z>-  - CorTNz,PNz  )_ 


< U)|l,N2  |‘  < w|z  \ l 


Z ( 


which  proves  the  claim.  By  the  Lumer- Ph i 1 1 i ps 


Iheorem  (cl.  ( 1 (>  ) ) .y ' is  the  infinitesimal  generator  of  a 


C'0  - semigroup  SN(t),  t >_  0.  Invariance  of 
_ N N NI 

a consequence  ol  .*  Z'  C Z . In  order  to  establish  (3.1)  it 
onl>  lemains  to  verily  that  the  hypotheses  (i)  and  (ii)  of  Lemma 
3.1  hold.  Condition  (i)  follows  from  our  comments  above.  In 
considering  condition  (ii),  we  choose  L/1  = Vk.  Bv  part  c)  of 
Lemma  2.2,  L/  and  (V-\I)i/k,  for  X sufficiently  large,  are 
dense  in  2 and  therefore  also  in  Z i . l or  J y we  have  the 

es  t ima  t e 


• M i ,N , 

under  S (t)  is 


M a . ' . iN  V \ W ' M «s  . 

q/  < | !»■  ^ rVv|..  ♦ |pVy-Vy|„ 

* k!  o £ 

g g J, 


< | iyi’  y iW  1 , ♦ |r‘^'sy'y 


(3.2) 


The  second  term  on  the  right-hand  side  of  (3.2)  approaches  zero 

as  N -*■  00  by  condition  (iii)  - (a)  of  our  hypothesis  above.  If 
N-  n N 

we  write  1’  41  = (4>  (0),4*  ) then  the  first  term  on  the  right-hand 
si ile  of  (3.2)  is  |(L(4>^),Dtl^)  - (L(ip),Di|i)|7  which  tends  to  zero 


.••v.k  .«i  . . *r. 


as  N 

fr  no 

by  ( i i i ) 

(b ) . Note,  t ha  t 

convergence  in  Z is 

equ  i va  1 ent 

to  convergence  in  2 . Therefore,  we  have 

1 i m j'' 

N *•■" 

NJ  - 

.i  for  a 1 1 

* k 

1 H - , 1 . e . 

hypothesis  (ii)  of  Lemma 

3.1  is 

sa  t i 

s f i cd . 

Roma  rk 

3.1. 

V 

In  order  to  prove  1.(4/  ) 

l.(T)  it  suffices  to  show 

• N * P 

i n 

C (with  sup 

norm) , because 

1.  is  a continuous  func- 

t i ona  1 

C * 

Rn  . 

Roma  rk 

3.2 . 

Since  y.  , 

k = 0,1.2 

i s dense  in  Z , it  is 

c 1 c a r 

that 

(iii)  - (a) 

ho  Ids,  if  1 i in  P‘ 

> = J for  all  Jr 

N 


for  some  integer  k ' 0. 

Remark  3.3.  The  proof  of  the  Trot  tor- K.ito  Theorem  also  yields 

w 

estimates  tor  the  rate  of  convergence  of  S (t)z  -*  S(t)z,  at 

least  tor  certain  z £ Z.  This  can  be  seen  as  follows.  Inequalities 

(4.2)  and  (4.3)  in  (10  , p.88]  imply 

| [SN(t)-S(t)]R(X0;.<y)2z|z 

< MetoT|  [R(\0;»  N) -R(X0  ley’’)  1 RO0  ; tv)  2 I z 
T 

♦ Mo"‘r  {j(R(X0;,VN)-R(X0;^)lS(o):|,do 

♦ ! tR(X0;  ^N)-R(XQ;  vy))S(t)R(X0;vi  ) = lz 


He  re 


for  t £ [ 0 , T 1 , T - 0, 


z t 


is  a fixed  real  number  such 


tnat  Aq  > io  and  R( A ; jy)  , R( A jj.*' ” ) is  the  usual  notation  for 
(XI-jy)  1 and  (XI-jylNI)  l,  respectively.  We  also  make  use  of 
the  estimates 


|S(t)  | 1 Me*1*,  |SN(t)  | < Mewt,  t>0 

and  |R(A0;.cyN)  | < ^ — with  some  constant  M _>  1 . liquation 
(4.14)  in  [in,  p . 9 1 J gives  the  estimate 

I [R(*0  ;jyN) -R(x0  ijy ) ] z 1 1 1 I [-^N-J^]R(X0;jy  )z  | z 

for  z E Z.  Thus,  for  some  constant  M = M(T)  we  finally  obtain 

| [SN(t)-S(t)]R(A0;jy)2z|z 

1 M { | [J^N-j^]R(X0;j^)2z|z 
T 

+ | I [-fl/N-J^]S(o)R(X0;_Qy)z|zdo 

+ I [j3fN-J3f]S(t)R(X0;j^)2z|z  j 


for  all  tE[0,T]  and  z E Z.  In  order  to  use  this  estimate 

2 

suppose  that  i1  is  a subset  of  ) such  that  for  positive 

constants  v and  p we  have 


18. 


lor  .ill  vi  4.  Lot  * | ( * bo  Miih  that  ; i 4^  implies 


.itul 


SCt ) 2 • 4 for  all  t t [0 ,T] 


S(t)  C \Q  1 - y)  z v 4 for  all  t t [0,T|. 


I lion  tor  each  ; » <r  j there  exists  a constant  v ■ v(:)  s 0 
such  that 


SNC t ) Sit)): 


* : N ‘ 


1.2 


for  all  t i [0,1]. 

In  the  case  of  spline  approximations  (.which  are  considered 

in  section  4),  it  is  not  difficult  to  find  subsets  4 and  4 ^ 

such  that  all  the  assumptions  given  above  hold  and  thus  convergence 

rate  estimates  as  mentioned  above  are  obtained.  lor  example, 

define  for  k s 2 the  sets 


y k - |i*(0),*)  i v-  k|i>,  + 1no)  = l i p 1 4>) , i o,i k 

I lien  (or  the  first  order  splines  discussed  in  Section  4 we  max 

choose  4 S'  a n d ■ y while  for  the  cubic  splines 

discussed  there  we  max  take  4 V*  and  4,  V ^ . However. 

I 

the  numerical  results  stronglv  indicate  that  the  estimates  one 
gets  this  wax  are  not  sharp  at  all.  Cert  a in  lv  this  problem 
needs  further  investigation. 


J 


m,t)  t 


19. 


\\c  next  turn  to  the  nonhomogeneous  problem  and  define  for 

- and  ft  I,,  . 

2,loe 


(t;n,<fr,f) 


,.N 


s'’(t)r’(M)  ♦ I s"(t-o)i,N(f  (o)  ,0)do, 
0 


(3.31 


t 2 °.  N = 1,2,... 


Theorem  3. 2 . 

a)  ho r (n.'fl  t Z and  f C I,,  j ^ we  have 

lim  =N(t;n,$,f)  = z(t;n,<fr,f),  t > 0,  (3.4) 

-*•<*>  — 

uniformly  far  t in  bounded  intervals. 

b)  For  any  t { > 0 the  limit  in  (3.4)  is  uniform  with 
respect  to  t C [ 0 , t j ] and  f in  hounded  subsets  of 
I-:(0,t1  ;r”)  . 

c)  I_f  (l  } is a sequence  in  I. , (0 , t j ; Rn)  converging 

weak  l v to  f , then 

lim  =N(t;n,d>.tk)  = z(t;n,4>,f) 

N , k-vt» 

uniform!  y for  t£[0,t.]. 


A 


Proof.  Of  course,  part  a)  is  a c 
Theorem  3.1  and  |SN(t):|„  < ewt 

N = . it  is  clear  that 


l im  sN'(t)pN(n.*)  = S(t)  (n,*) 

^ VO© 


onsequence  of  part  b).  from 

I : I ■*  for  t > 0.  Z ( 1 an. 


on  i t o i ml v for  t t [ 0 , t j J . 
Following  [2  I Je 


•rN,  . > „n 


T (t)  : R11 

- Z bv 

TN  ( t ) r. 

S 6 R , 

t > 0. 

Then  it 

for  any 

t > 0, 

M 

T (t)  -> 

for  anv  t > 0.  T»(t)  - T(.)  as  N . - in  the  uniform  operator 
norm  anti  j V(o)-T(o)  |-\lo  - 0 as  N - Thus  for  t6|0.t.| 


[s  (t-o)P‘  (f Co) ,0)  - S(t-o) (f (a) ,0) ]do 


i I'1'  (t-a)  - T(t-o)  | | f ( a ) | do 


(0)  - T^i-’-r,->itlo.,1;R") 


which  proves  part  b)  of  the  theorem. 

Part  c)  follows  from  the  inequality 


z‘  ( t , fk)  - z (t  ,n,$,f ) L 

M 

< I 7 ^ f t • n r 


- ,z  Ct;n’*’f  ) ‘ zCt;n,d',fk)  | z ♦ |z(t;n,*,fk)  - z(t;n,$,f 


k 

using  part  b)  above  and  the  fact  that  {f  } is  in  a bounded 

n 

subset  of  L2(0,t^;R  ),  and  b)  of  Lemma  2.4. 


Remark  3.4.  Theorem  3.2,  part  c) , is  of  special  importance  for 
the  application  of  our  approximation  results  to  optimal  control 
problems.  See  the  discussions  in  [2]  and  [3]. 

Remark  5.5.  If  we  write  zN  ( t ; n , 4> , f ) = (xN(t)  ,yN(t) ) , t >_  0 , 
then  (3.4)  implies 

lim  xK(t)  = x(t;n,<J>»f) 

N-*-oo 

uniformly  for  tetO.t^,  > 0.  This  is  clear  by  definition 

of  the  norm  in  Z. 


In  order  to  obtain  algorithms  which  can  be  implemented  on  a 
computer  we  assume  from  now  on  that 


,N 


dim  Z = < 00 » N = 1,2,... 


Then  zN  (t ; n , <P , f ) given  by  (3.3)  is  the  unique  solution  of  the 
ordinary  differential  equation 


. M W M W 

z (t)  = .f/  z (t)  + P1N(f(t)  ,0)  , t > 0 


(3.5) 


_N 

in  Z . 


w N 

z (0)  = PN(n,<P) 


We  fix  a basis  for  . Since  C 


:N 


,N 


,N 


N 


,..1.2 


\%e  have  6.  = (B'.  (0)  , B'- ) , j = 1 kM,  with  B • 6I(  ’ . Define 


J 


J 


the  n * k.,  matrix  function  6 

N 


by 


MM  M 

s = cb”,  — 


N 


and 


-'N  N M 

B = (6  (0) ,B  ) 


N N 

Then  any  z £ Z can  be  written  as 


N sN  N f J,n.  N CN  N, 
z = B a = (8  ( 0 ) u , 3 01  ) 


N N N ^ m 

where  a1  = col  (a‘,  , . . . , a‘,  ) £ R is  the  coordinate  vector  of 

1 kN 

N 

z with  respect  to  the  chosen  basis.  The  matrix  representation 
N . N 

of  -0/  restricted  to  Z with  respect  to  this  basis  is  denoted 

by  AN,  while  w^(t)  and  F^(t)  = col (F^(t) , . . . , F^  (t))  are 

N 

the  coordinate  vectors  of  the  solution  z^(t)  of  (3.5)  and 
PN(f  (t) ,0)  respectively.  That  is,  z^(t)  = B^w^(t)  and 
P^(f(t),0)  = B^F^(t).  Then  system  (3.5)  is  equivalent  to  the 
system 

wN(t)  = ANwN(t)  + FN(t),  t > 0 


N N 

w (0)  = wQ 


(3.6) 


M W M MNJM 

in  R Where  8 wQ  = P ( n , J> ) . If  w ( t ; wj  , F ) denotes  the 
solution  of  [3.6),  then  by  Theorem  3.2  we  have 


lim  BN^N'(t;wp,rN)  = (x  (t  ;n,4>,  f) , x.(n,$,f)) 

^-►oo  L 


uniformly  on  bounded  t- intervals  and  uniformly  with  respect  to 
f in  bounded  subsets  of  L, (0 , t j ;Rn) , ^ > 0.  According  to 

Remark  3.5  we  also  have 


lim  8N(0)wN  (t  ;w^,FN)  = x(t;n,<J>,f) 

uniformly  in  t GlO.tj],  tj  > 0. 

In  order  to  solve  system  (3.6)  on  a computer,  one  must,  of 
course,  know  how  to  compute  PN(n,<(>)  and  the  matrix  AN.  We 
first  show  how  to  compute  the  coordinate  vector  of  PN'(  n , <#>)  for 

(n,tf>)  £ Z.  Since  P'N'  is  the  orthogonal  projection  Z^  -*•  ZN  tlie 

element  p\n,4>)  is  uniquely  determined  by  the  orthogonality 

relationship  (in  Z^) 

{PN’(n,<»  - (n,0))  1 ZN. 


This  is  equivalent  to 

-'MM 

<b  « - (n  ,'f)  ,bn)z  = 

g 


o 


24. 


or 


„N  N 
Q a 


hN(n,4>) 


(3.7) 


where 


QN=(6N,6N)z  = BN(0)T8N(0)  + f BN(0)T6N(0)g(0)d0 

g ‘ -r 

and 

hN(n,4>)  = <6N,(n,4>))7  = BN(0)1n  + [ BN(0)  r<}>(0)g(0)d0. 

g 1 -r 

Thus,  we  see  that  in  order  to  get  the  coordinate  vector  of 

*r'N (n , 4>)  we  have  to  solve  (3.7).  Note  that  (Q^j  k must  exist 
since  P^(n,<j>)  is  uniquely  determined  by  C n , ) - 

With  respect  to  F^(t)  we  have  h^(f(t),0)  = 8^(0) ^f (t) 
and  therefore 


QNFN  ( t)  = BN(0)Tf(t). 


Finally,  the  matrix  is  calculated  in  the  following  manner. 

r IN  .Ns  c 7N  , N c DkN  . N _ ,,kN  . 

For  <p  = (<P  ( 0 ) , <t>  ) E Z define  a E R and  y E R by 


A 


N 


BNaN 


and 


jy  <t> 


§nyn 


Now,  _ayN$N  = P Njy$N  = PN(L(4)N)  , D4>N)  . In  light  of  the  calcula- 

N 

tions  given  above  we  see  that  y is  the  solution  of 


I 


>3 


... 


25. 


J N , N,,  , .N-.  n.N. 

Q Y = h CL(<|»  ),  D<*>  ) . 


N, 


But  by  linearity  of  the  map  (n,4>)  h (n,4>)  and  L we  have 
tha  t 


,Nri  r.N,  n.N,  ,Nn  , J N ,noNs  N-. 
h (L(<f>  J , Dip  ) = h (L(B  a ) , (D$  ) a ) 


,,N  N 
= H a 


where 


N N N N 
llN  = hN(LC3  ),DBN) 


= BN(0)TL(BN)  + (°  6N(0)T(DBN) (0)g(0)d0 

■*  - r 


(3.8) 


Thus  , we  obtain 


AN  = (QN)'1HN. 


Of  course,  on  a computer  one  never  actually  computes  (Q^)  * but 
rather  solves 


„N  N 

Q Y 


„N  N 
H a 


(3.9) 


V « ^ 


directly  in  order  to  obtain  = A^a^. 


26, 


4 . Spline  approximations 

In  this  section  we  show  that  the  general  scheme  given  in 

N 

Section  3 can  be  realized  by  choosing  the  Z as  certain  sub- 
spaces of  spline  functions.  In  order  to  keep  the  presentation 
brief  we  give  the  details  only  for  first  order  splines. 


.N 


Corresponding  to  the  partition  = -jjjj,  j = 0,...,N, 

of  [ - r , 0 ] we  define  Z^1  = j (4>  (0) , <}>)  G &\<p  is  a first  order  spline 

function  with  knots  at  t j , j = 0,...,n|.  A first  order  spline 

function  <f>  on  [-r,0]  with  knots  at  the  {tj}  is  simply  a 

continuous  function  [-r,0]  which  is  linear  on  each  subinterval 
N N N 

[tj’tj-l^*  j = Let  P^  be  the  orthogonal  projection 

Z - Z?  and  = P^-crf,  N = 1,2,... 

g 1 111 


Theorem  4.1.  The  approximation  scheme  {Z^,P^,_of  satisfies 
all  conditions  of  Theorem  3.1  and  dim  Z^  = n(N+l). 


Proof . Conditions  (i)  and  (ii)  of  Theorem  3.1  are  trivially 

N 

satisfied.  It  is  also  clear  that  dim  Z^  = n(N+l).  For  the 

a 2 

verification  of  condition  (iii)  we  take  k = 2.  We  fix  ip  G 

M ~ M 

and  put  = P^.  Let  then  denote  the  interpolating  first 

order  spline  function  defined  by 


^(tj)  = ’Ktj) , j = o 


9 


N. 


Using  the  well  known  convergence  properties  of  interpolat- 


27. 


in£  splines  (cf.  for  instance  [11,  Thm.  2.5])  and  the  fact  that 

i$-$i  wc  obtain  immediately  that 
TS 

l^-^lz  1 I I z 5 | V ^ < 0 ( '-v)  (4.1) 

K g <•  ' N“ 

as  N -*■  ®.  Thus,  condition  (iii)  - (a)  of  Theorem  3.1  is  satis- 
fied (cf.  also  Remark  3.2).  The  theorem  in  [ 14  J quoted  above 
also  provides  the  estimate 


|D(^-4>)  l2  i 0(g). 


Therefore 

1 1U^N-^)  1 2 1 |D(^-*)|2  - |1)(^-U»N)  l2 

(4.2) 

5 o(g)  - 

where  for  the  second  term  on  the  right-hand  side  we  have  used  the 
Schmidt  inequality  for  polynomials  of  degree  one  on  each  of  the 
intervals  [tj.tj.j_J,  j = 1 N (cf.  [M  , Thm. 1.5]).  IVe  ob- 

serve that  the  constant  M in  this  estimate  is  not  dependent  on 
N.  Using 


i * kN-q,  < o(ij) 


28  . 


(note  that 
f i na 1 1 y t he 


cst imate 


for 


4 C C) 


along  with  (4.2),  we  obtain 


|IH*N-H»)|2  < 0(1) 


(4.3) 


as  N -*•  <».  This  proves  the  second  requirement  in  condition  (iii)  - (b) 
of  Theorem  3.1  and  it  only  remains  to  verify  that  l.(TN)  -*•  I.(tj') 
as  N ■*■«>. 

For  0 £ (-r,01  we  have 


VN(0)  = N ( 0 ) 


r n 

(PV  ) (o)do 

0 


and  hence 


|fNl°)  * * (°) I 5 |fN(0) 

< |tN(0)  - no) | 

< lAo)  - ♦ (0)  I 


- V ) (o) ! do 


V ) (o) | "do 


1/2 


This  estimate  together  with  (4.1)  and  (4.3)  implies 


|A0)  - v (Q)  | s O(^) 

as  N -*•  uniformly  for  0 o l-r,0).  Therefore,  |L(Y|N)  - L (V ) | < 
O(^j-)  (see  Remark  3.1)  and  the  proof  of  Theorem  4.1  is  completed. 


Remark.  4.1.  The  estimates  given  above  reveal  t ha 


N - 

& ^ 


as  N 


08  lor  ^ ( '*'■ " • Based  on  this  estimate  the  rt 


1 1 e o f 


convergence  tor  S (t)  indicated  in  Remark  3.3  is 


Is  It  u'  S(  t if  I : . 0(A) 


a s N -»  <« , Co  r if  i n V ' . 


Adhering  to  the  general  outline  given  in  Section  3 for 

representation  ot  system  (3.S)  in  vector  matrix  form,  we  used  the 

following  coordinate  representation  for  our  numerical  calculations 

with  first  order  splines  reported  in  Section  5 below. 

N 

1 ^ *■  j • 3 B » • • • * N , denote  the'  scalar  first  o r d e r spline 

function  on  [ - r , 0 1 characterized  by 


N N 

cj  ^ j)  = 6ij  • 3 = 0 N, 


6ij  being  the  Kronecker  symbol.  Then  the  matrix  pN  is  given  by 


,N  N N 

* = o « 1 


where  W denotes  the  Kronecker  product  and  1 is  the  n x n 
identity  matrix.  An  element  z in  Z,  with  coordinate  vector  oN 


where  the  vectors  E Rn  are  such  that  = colCaJJ , . . . ,aJJ)  . 
For  the  case  where  m = 1 in  (2.2)  (where  we  have  g(0)  = 1), 


simple  calculations  yield 


n 

i 


0 N , 


j » ( A(0)Cj  (0)d0  , j = 

' r N 

max  ( t , , - r I 
J + 1 


and 


for  N = 2,3,...  . 

N 

Finally  F [tl  is  given  by 

N N - 1 N T N 1 

F (t)  = (Qj)  P (0)  f c t ) = CQj)  collfftl,  ..,01. 

It  is  the  happy  circumstance  that  all  the  ideas  detailed 

above  for  first  order  splines  carry  over  to  splines  of  higher  order 

For  instance,  in  the  case  of  cubic  splines  one  may  take  k = -1  in 

order  to  prove  a theorem  analogous  to  Theorem  4.1.  Instead  of 

Theorem  2.3  in  (14]  we  have  to  use  Theorem  e.7  in  | 1 4 ] . 

N 

It  can  be  seen  that  dim  I,  = nfN+31.  For  the  construction  of  a 

.> 

N 

basis  of  Z~  we  take  the  cubic  de  Boor  splines  (see  [14,  p.  73  ]1. 

N M 

The  matrix  Q,  is  a 7 band  matrix,  H, 

1 s 


is  a P band  matrix,  which 


3 2 . 


we  shall  not  i;ivc  here  explicitly  (but  see  [13]  for  the  details), 
lor  cubic  splines  we  obtain  the  estimate 

| - a/y  U < O(i^) 

L ' N 

~ .4 

as  N -*■  00  for  i>  6 . Use  of  the  estimates  in  Remark  3.3  will 

y i e 1 d 

| ( t ) v - s ( t ) v | < o(-L.) 

N' 

as  N -*■  00  for  Y in  S/  \ 

5 . Examples 

In  this  section,  we  discuss  some  of  the  examples  we  used 
in  computations  in  order  to  investigate  typical  features  of  the 
spline  approximations.  For  comparison,  we  also  give  the  results  for 
the  so-called  averaging  approximations  which  are  discussed  in 
( 2 ] and  [ 5 ].  In  the  following,  we  write  AY , Sj  and  S, 

for  the  averaging,  first  order  and  cubic  spline  approximations, 
respectively.  Examples  1,2  and  5 for  AV  and  Sj  were  computed 
on  the  HIM  360/67  at  Brown  University,  whereas  Examples  3,4  and  3 

for  S_  were  run  on  the  IINIVAC  1100  of  the  Rechenzent rum  Graz.  Kc 

¥ 

are  grateful  to  Dr . D.  Reber  and  Mr.  P.  Rotter  for  their  assistance 
with  the  calculations.  Our  main  interest  in  performing  the 
calculations  was  to  demonstrate  that  the  algorithm  presented  in 
Section  4 is  numerically  feasible  and  to  obtain  information  about 


"*  — — - .ccur.cy  of 

'"'OSM"S  « — «re  „„t  opti<lj;i.1  ;P|PrOXiM,io"-  Therefore.  ,„o 

S 0raK"'  °tC-  routue,  f„  ' rCS,'l’Ct  10  «■*«<»«  ti,„, 

h^dst s Iu;;;o";n*  ^>.  ,,6) 

'"d°r  eI«ination  and  fourt|] 


Table,  below  , 

is  »»«  Cappro,i.„etv  cho  " "'C  f°U°"'i"s  »«tati0„.  s« 



‘"•■'t  segment  of  j , ... „ * where  WN  ( t ) 

in ons  1 on  n in  e.  r ’ 15  the 

C°1  ( t ) N f he  -s° Inti  on  voc  t , V 

1 Vim)  of  the  approY.  ctor  w\t)  - 

,"Stan"  U.  p.  ;'""S  °"T  (see  f0r 

' * ' 1 n case  nf  w 

^ o 1 u ^ » on  o • » . • ' » 3n tl  x / 1 ) 

• Slm“arly,  SN  of  „ '(<)  -s  the  true 

tho  differences  1 S3  “ ,h?  appropriate  n.,.  , 


i-  or  ,V’ 

S''  S thC  ^Topria.e  chojce  fron 


/(f«V(0)wNl  Ct;wN  ,:N 

J “•  1 ■ xj<':n,«,ni  , 


J = 1 


Crrors  for  the  s or  $ 

In  this  examoie 

osci"ator  ^ ^ r..;or7n";.T ,hc  equation  f°r  ■ --ped 

8 ^ 3nd  «“*■«  forc, 


x(t)  + x(t)  ♦ 


Wlth  initial  c 


om!  it  ions 


•Mt-l)  = 10 


x(0) 


COS  V 


* *(#)  = 


-l  < 0 < 


The  solution  on  the  interval  (0,2) 


is  given  by 


x(t)  = u ( t ) -9  - sin  1 ♦ lOt  ♦ (10  + \ sin  1 - 4 cos  lie'1 

♦ !,  (sin  l - cos  l)sin  t ♦ 4 (sin  1 + cos  l)  cos  t 

for  t E ( 0 , 1 1 and 

x(t)  = u(t)  - 29  - 2 sin  1 ♦ cos  1 ♦ (19  ♦ sin  l)(t-l)  - 5 C t - 1 ) 2 

♦ (4-  * T Sin  1 - cos  l)e'(t'n  ♦ 

♦ (10  ♦ 4 sin  1 - 4 cos  l)(t - 1 ) e “ (t * 1 1 

♦ 4 (sin  1 • l)cos(t-l)  + 4 (1  - cos  l)sin(t-l) 
for  t E (1,21  . 

For  this  example,  t he  calculations  were  carried  out  for  AY 
and  Sj.  Tables  1 and  2 show  the  numerical  results  for  x(t)  and 
x(t).  As  is  the  usual  practice,  this  second  order  equation  was  con- 
verted for  computational  purposes  into  a 2^2  system  of 
first  order  equations  for  the  unknown  functions  x ^ ( t ) = x(t)  and 

x,(t)  = x(t).  The  results  clearly  show  that  the  convergence  for 
1 - £ 

AV  is  like  1/N  , with  c a small  positive  number,  whereas  for 

Sj  it  is  like  1/N".  Another  typical  feature  exhibited  by  this 
example  is  that  the  relative  error  for  AY  is  increasing  much 


faster  with  time  (in  this  example  for  AV 


the  relative  error  for 


t 


vs  about  seven  times  the  relative  error  at  t 


it  i s for  S . 
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0.25)  t ban 


luble  1 


r~  r 

t 

— r 

x(t) 

4 

6u; 

:W 

4 

4 

4 

T 

! 0.25 

1 .270-18 

0.0012b 

0 . 00 0b 4 

0.00032 

0.005O4 

0.00132 

0.00033 

| 0.5 

1 . 993b" 

0.0037" 

0.00194 

0.00098 

0.00895 

0.0022b 

0.0005b 

I °-"S 

5. Ob 148 

O.OObOl 

0.00300 

0.00151 

0.01170 

0.00292 

0.000" 3 

; i-o 

4.39272 

0.00951 

0.00421 

0.00189 

0.0134b 

0.00337 

0.00084 

1.25 

5.92593 

0.0226b 

0.01089 

0.0052" 

0.014"8 

0.00374 

0 . 00094 

1.5 

".60007 

0.05099 

0.02b  lb 

0.0135* 

0.01281 

0.00307 

0.00075 

1.75 

9.34402 

0 . 09 1 1 7 

0.04813 

0.02484 

0.00"b3 

0.00200 

0.00050 

2.0 

— 

11.08 330 

0.1375 



0.0725 

0.0373 

0.0043 

0.0010 

0.00030 

Table  2 


t 

un 

4 

8 1,1 

AY 

_ 

4 

_jL_ 

4 

4 

0.25 

2 . 0b9b9 

0.00882 

0.00453 

0.00230 

0.00" 21 

0.00221 

0.00050 

0.5 

I 

3.64428 

0.01006 

0.00521 

0.002b" 

0.0052b 

0.00120 

0.00028 

! 0.75 

4.84445 

0.00851 

0.00319 

0.00135 

0.00156 

0.00057 

0.00014 

|,.o 

5.76581 

0.02570 

0.01038 

0.003S9 

0.00083 

0.00011 

0.00002 

1.25 

6.45956 

0.0832b 

0.04532 

0.02428 

0.00054 

0.00003 

0.0000" 

1.5 

6.88559 

0.1406b 

0.0"b9S 

0.04057 

0.009" 7 

0.00319 

0.00063 

1.75 

7.01599 

0. 1"b85 

0.09403 

0.04835 

0.00853 

0.00184 

0.00045 

2.0 

6. 849"  2 

0.19072 

0 . 09909 

0.05029 

0.00b 11 

0.00160 

0.00044 

3o. 


fix  a nip  1 e 2.  Here  we  deal  with  the  equation 
delayed  damping 


for  t he  oscillator  w i t h 


xit)  + x(t- 1)  + x(t)  = 1 


with  initial  da t a 


x = *(0)  = o for  0 t [-1,0]. 

The  solution  on  [0,2]  is 

x(t)  = 1 - COS  t for  t C [0,1] 

and 

= * ‘ cos  t + 3 (t  - 1 ] cos  ( t - 1)  - 4 siiut-n  for  t E [1,2]. 

Again,  t tie  calculations  were  done  for  AY  and  Sj.  The  data  given 
in  Tables  3 and  1 depict  behavior  similar  to  that  found  in 
Example  1. 


i .1 M o 


0.03  H'S8 


0.12:il ' 1 0.000016 


0 

0.000661 

0.00015' 

0 

0.001185 

0.00029S 

0.000010 

0.001609 

0.000401 

0.000580 

0.001767 

0.000444 

0.003529 

0.001985 

0.000481 

0 . 0068 1 ~ 

0.001223 

0.000293 

0.009290 

0.000260 

0.000070 

0.010383 

0.00188' 

0.0004 '7 

Table  4 


0.247404 

0.000001 

0.479426 

0.0002S 3 

0.681639  j 

0.004566 

0.841471 

0.024575 

0.918059 

0.044008 

0.877639 

0.043160 

0.728371 

0.028576 

0.488562 

0.001830 

0.000017  0.000001 

j 

0.001231  ’ 0.000214 

0.013283  0.006090 

0.025291  0.013826 

0.023308  0.011 S91 
0.014968  0.007548 

0.000834  0.000551 


0.001138  0.000314  0.000074 
0.000913  ; 0.000218  0.000054 
0.000148  0.000081  0.000018 
0.000642  0.000216  0.0000'S 
0.001180  0.000122  0.000041 
0.002876  0.000S54  0.000202 
0.004395  0.001004  0.000255 
0.004057  0.001039  0.000265 


Example  3.  We  consider  a 2-dimensional  system  described 
by  the  equation 


38. 


x(t)  = 


0 

1 

O 

O 

x(t)  + 

-1 

0. 

.0  -l. 

x(t-l) 


with  initial  condition 


x(*>)  = col (0 , s in  2 nO ) for  0 E [1,0]. 

The  solution  on  [0,2]  is  given  by 

xi  (t)  = — ( - 2 it  sin  t + sin  2nt), 

(2*)  -1 

x2(t)  = x1(t) 
for  t G [0,1]  and 

XjCt)  = — ~~2 — (*  sin  1 * sin(t-l)  - — cos(t-l) 

C 2 tt  ) - 1 < (2tt)2-1 

+ j (t-l)sin(t-l)  + — cos  2ntJ, 

x 2 ( t ) = x x ( t ) 

for  t E [1,2].  In  Tables  5 and  6 we  give  the  numerical  results 
for  AV  and  S^.  For  AV  we  have  convergence  like  1/N^  C, 
o > 0 small.  The  data  for  cubic  splines  do  not  behave  as  regularly 
as  the  data  for  AV  or  as  that  for  in  the  previous  examples. 

A possible  explanation  for  this  is  the  more  complicated  structure 
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of  the  mat r ices  which  one  encounters  in  loi  instance,  one 

must  he  very  careful  when  deal i hr  with  equation  (3.9).  We  note 
t|iat  N = 4 for  S , Rives  results  comparable  to  that  for  N = 12S 

for  AV . 


Example  4.  We  next  use  an  example  (clue  to  l’opov)  for  a degenerate 
system  where  we  have  (1  ,-2,-1)  x(t  ; n,*^)  = 0 for  t ' 3 and  all 
(n,4>)  G Z.  The  equation  is 


x(t ) = 


0 

p 

(1 

0 

0 

0 

0 

0 

- 1 

X (t)  + 

1 

0 

0 

0 

0 

0 

0 

-y 

(1 

x(t  - 1) 


We  choose  the  discontinuous  initial  data 

n = col  (1,1,1)  and  4>(°)  = 0 for  Of  [1,0) 

The  solution  is  given  by 

x(t)  = col(l+2t-t~,l-t,l)  for  t G [0,11, 
x(t)  = c o 1 ( 2 , 0 , 2 + 4 t - 1 “ ) for  t G [1,21 

a n d 


x(t)  = co  1(2, 0,2)  for  t > 2. 


N 

Here  the  initial  data  are  in  the  subspaces  Z lor 
for  S,.  Thus  we  in  ip, lit  expect  that  AV  will  initially 
for  t small)  give  a better  approximation  than  S,. 
simple  algebraic  structure  of  the  approximating  ODD ' s 
AV  it  is  not  difficult  to  prove  that  we  have 


AV  but  not 
C i -e. 

Due  to  the 
in  case  of 


1 im  w^(t)  = col (2  + 1 , 0, 

£ >co  l N 


for  all  N = 1,2 So  for  AV  we  also  can  expect  a very  good 

approximation  of  x_(t)  for  t large.  The  error  for  x.(t) 
and  x,(t)  should  tend  to  ~ as  t 00 . That  this  is  true  is 
shown  by  the  numerical  data  displayed  in  Tables  7,  8 and  9.  Again, 
a characteristic  feature  of  S.  is  that  the  error  is  not  in- 
creasing in  magnitude  over  the  time  interval  under  consideration. 
Note  that  S , for  N = 4 gives  a better  approximation  of  x^(t) 


for 

t > 1 and  of  x ^ C t ) 

for 

t > 2.6 

t ha  n 

AV 

for  N 

= 128. 

With 

respect  to  x,,(t)  we 

can 

sec  that 

S3 

for 

N = 4 

gives  a 

better  approximation  than 

AV 

for  N = 

128 

around  t = 

1 where 

the  true  solution  x,(t)  has  a jump  discontinuity  in  the 
derivative.  That  AV  behaves  better  for  t < 1 - e and 
t > 1 + e,  e some  positive  number, is  due  to  the  specific  choice 
of  the  initial  data. 


. 
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Table  7 


t 

- - - i 

XjU) 

r32 

6av 

*64 

6av 

.128 

6AV 

4 

3 

& 

4* 

1 

0.2, 

1 

36000 

0 

0 

0 

0.01211 

0.00493 

0.00247 

0.4 

j 

1 

64000 

0.00001 

0 

0 

0.00885 

0.00528 

0.00001 

0.6 

1 

84000 

0.00001 

0 

- 

0 

0.00446 

0.00334 

0.00190 

0.8 

1 

96000 

0.00105 

0.00015 

0.00001 

0 . 00980 

0.00074 

0.00146 

1.0 

2 

0 

0.01416 

0.00729 

0.00372 

0.00020 

0.00729 

0.00372 

1.2 

o.o: 

866 

0.01519 

0.00777 

0.00725 

0.00424 

0.00175 

1.4 

0.03100 

0.01562 

0.00781 

0.00160 

0.00460 

0.00130 

1.6 

0.03123 

0.01563 

0.00781 

0.00527 

0.00406 

0.00085 

1.8 

0.03 

125 

0.01. 

563 

0.00 

781 

0.00331 

0.00281 

0.00242 

2.0 

0.00256 

0.00117 

0.00130 

2 . 2 

0.00306 

0.00049 

0.00096 

2.4 

0.00053 

0.00178 

0.00214 

2.6 

0.00212 

0.00246 

0.00115 

2.8 

0.00054 

0.00251 

0.00092 

3.0 

2 . 

0 

0 . 03 

125 

0.01563 

0.00781 

0.00117 

0.00200 

0.00191 
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Table  8 


t 

X2 

(t) 

<5 

32 

AV 

6 

64 

AV 

6 

128 

AV 

,4 

6s 

53 

*8 

S3 

.16 

6S 

b3 

0.2 

0 

.8 

0 

0 

0 

0.00947 

0 . 00536 

0.00230 

0.4 

0 

.6 

0.00001 

0 

0 

0.01206 

0.00452 

0.00035 

0.6 

0 

.4 

0.00017 

0 

0 

0.00526 

0.00260 

0.00202 

0.8 

0 

.2 

0.00881 

0.00200 

0.00020 

0.00153 

0.00015 

0.00139 

1.0 

0 

0.07034 

0.04980 

0.03524 

0.02133 

0.01895 

0.00949 

1.2 

0 

0.01375 

0.00384 

0.00056 

0.01505 

0.00168 

0.00001 

1.4 

0.00159 

0.00009 

0 

0.01369 

0.00068 

0.00094 

1.6 

0.00012 

0 

0 

0.00105 

0.00078 

0.00208 

1.8 

0.00001 

3 

0.01148 

0.00105 

0.00193 

2.0 

0 

0.00444 

0.00123 

0.00020 

2.2 

0.00635 

0.00108 

0.00173 

2.4 

0.00521 

0.00066 

0.00202 

2.6 

0.00205 

0.00019 

0.00073 

2.8 

0.00407 

0.00021 

0.00084 

3.0 

) 

) 

0 

0.00039 

0.00051 

0.00119 

Table  9 


x.tt) 


0.01131  0.00517  0.00235 


0.01278  0.00171 


0.0002. 


0.0051 


0.00028  0.00373  0.00202 


0.01 058  0.00386  0.00030  0.00-172  0.00100  O.OOlo-J 

I 

.0  0.12052  1 0.09231  0.06676  0.04871  0.03399  0.01807 

| 

.30  0.00115  0.00751  0. 00005  0.02598  0.00193  0. 00149 

.01  0.02780  j 0.01514  0.00781  0.022OO  0.00337  0.00228 

,84  0.03040  0.01559  0.00781  0.00015  0.00168  0.0024S 

.96  0.02546  0.01424  0.00760  0.01932  0.00086  0.00153 

.0  0.00207  0.00073  0.00026  0.01007  0.00266  0.00049 

0.02246  0.01347  0.00746  0.00901  0.00349  0.00223 


0.02950  0.01547  0.007S1 


0.03098  0.01562 


0.03122  0.01562 


0.00375  0.00198 


0.00071  0.00321  | 0.00068 
0.00714  0.00211  I 0.00049 


.0  0.03125  0.01562  0.00'Sl  0.00465  0.00076  0.00059 


■lb. 


Hxamp 1 c 3.  For  t ho  scalar  equation 

x 1 1 ) = Sx(t)  ♦ xtt-1) 

with  initial  da  t a 

x (d ) = 5 fo r 0 t [-1,0] 

wo  Rive  the  numerical  data  for  AY,  8^  and  S.  in  fables  10  and  11. 
The  true  solution  is 

x ( t 1 = oo ^ 1 - 1 for  t 6 10,1] 

and 

xft)  = lx ( 1 ) - | ♦ bit- 1)}  eSU"  11  ♦ ~ 

for  t 6 [1,2].  The  initial  data  are  in  z\  n = 1,2,...,  for  all 

three  approximations  AY,S  ,S^.  The  numerical  results  show  that  the 
approximation  by  AY  is  initially  better  than  by  8^  and  8,.  This 
appears  to  be  the  case  in  all  examples  where  the  initial  data  are  in 


the  subspaces 

-N  .. 

Z for 

AY,  Sj  and 

s, . 

.•> 

But  the  relat 

ive  error 

increases  much 

fast er 

with  increasing 

t ime 

for  AY  than  for 

S1 

and  S,.  For 
.> 

i ns  t a nee 

, f o r N = 3 2 , 

S1  W 

ill  do  better 

t han 

AY 

shortly  after 

t = 2 . 

The  numerical 

dat  a 

a 1 so  ind i e a t o 

accumulation  of  errors 

for  S,  if  N 
3 

and 

t are  large 

. The 

r a t e 

of  convergence 

f o r S . 

again  is  1 i k e 

1/N 

For  S, 

3 

i t seems 

to  bo  like  1 / ! 

> 0. 

0 

0.0001 
0.0001 
0.011 
0 . 530 

•i.  r 

23.4 

108.7 

474. 


0 

0.2853 

1.5512 

7.352 

32.380 

130.14 

555.8 
2210.9 
8721  . 


0 

0.0 7 So 


1.937 
8 . 5 3 1 
3 5 . 9 3 
14  7.0 
588. 3 
2316. 


0 

0.0102 


0.4  000  0.105' 


0.402 

2.163 

9.11 

37.3 

149.4 

588. 


lab  1 e 11 


t 

xtt) 

0.2 

15.309601 

0.4 

43.334337 

0.6 

119.513222 

0.8 

326.588900 

1 . 0 

889.478955 

1 . 2 

2420.772761 

1.4 

6588 . S65S18 

1 .6 

17954 .15321 

1 . 8 

48815.25691 

2.0 

132871 . 3779 

0.030796  | 0.004862  0.000289  0.000450  0.000032 

0.051843  0.003000  0.001761  0.000128  0.000303 

0.260502  0.012450  0.001268  0.000888  0.001223 

0.781007  0.018429  0.000120  0.004050  0.0031S0 

2.388011  0.094175  0.002117  0.013687  0.008004  1 

7.268330  0.243622  0.011288  0.038706  0.076608  j 

I 

21.936935  0.700860  0.028040  0.155481  0.206237 

65.33753  2.00746  0.07668  0.54728  1.26327 

193.16251  5.69406  0.52984  1.80974  4.94164 

566.4524  16.8275  1.4527  4.0407  17.0516 


4 8 . 


From  the  numerical  results  presented  here  we  can  draw  the 
following  conclusions: 

(a)  The  rate  of  convergence  for  Sj  in  all  cases  was  like 

■> 

1/N  compared  with  1/N  or  less  for  A\  . Ihis  superiority  is  not  only 
true  for  examples  where  we  considered  just  an  initial  value  problem  hut 

also  for  examples  (.not  presented  in  this  paper)  where  optimal  control 
problems  and  identification  problems  were  considered.  For  S,  con- 
clusions to  be  made  from  the  data  are  not  as  clear.  In  case  of  the 
scalar  equation  of  Example  S the  rate  of  convergence  seems  to  be  better 
than  1/N*.  For  the  3-dimensional  system  of  Example  4,  N = lb  gives 
just  a slightly  better  approximation  than  N = 4.  A possible  explana- 
tion for  this  behavior  is  that  equation  (.3.9)  must  be  dealt  with  more 
carefully  than  we  did  in  our  preliminary  calculations. 

(b)  For  the  same  value  of  N,  the  accuracy  of  approximation  over  large 

t-intervals  for  S.  and  S„  is  considerable  better  than  for  AY . 

1 a 

The  latter  scheme  does  better  initially  if  the  initial  data  are 

N 

already  in  the  subspaces  Z for  AV  and  in  some  exceptional  cases 
where  the  simple  algebraic  structure  of  AY  is  advantageous 
(see  Example  4). 

(c)  For  Sj  implementation  of  the  algorithms  on  a computer  is 
almost  as  easy  as  for  AV.  For  S.  the  matrices  appearing  in  the 
algorithm  arc  more  complicated.  On  t he  other  hand  in  many  cases  the 
results  using  AV  were  comparable  in  accuracy  to  those  obtained  using  S. 

= 4 only  if  N > 100  was  taken  in  the  AV  approximation. 


with  N 
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